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PREFACE 


In recent years much attention has been devoted to the development of simple and computationally fast analyti- 
cal approximations to the radiative transfer equation. This has largely been the result of the need to parameterize the 
radiative properties of clouds and aerosols in general circulation climate models. In these and other climate model 
applications, it becomes necessary to rapidly calculate the plane albedo, total transmission and fractional absorption 
as a function of optical thickness and solar zenith angle for a wide range of atmospheric conditions. 

This report contains a comparison of the absolute and relative accuracy of eight different radiative transfer 
approximations as a function of optical thickness (0.1 <r t <100) and cosine of the solar zenith angle (0 </a Q < 1). 
Contour plots of the approximate plane albedo, total transmission and fractional absorption are presented for each 
model, as well as contour plots for the relative and absolute errors in each model. These results have been obtained 
for four values of the single scattering albedo (viz.., 1 .0, 0.99, 0.9 and 0.8) and for a cloud phase function having an 
asymmetry factor g = 0.843. The radiative transfer approximations considered in this report are asymptotic theory 
for thick layers and the following widely used two-stream approximations: Coakley-Chylek’s models 1 and 2, Meador- 
Weaver, Eddington, delta-Eddington, PIFM and delta-discrete ordinates. The baseline computations for these com- 
parisons were obtained using the doubling method. 

The authors are grateful to Howard G. Meyer for aid in performing the computations and in generating the 
numerous graphical results. The research was in part supported by NASA Grant NAG 5-309. 


iii 




CONTENTS 


Preface 111 

Summary ^ 

Introduction * 

Multiple Scattering Computations 2 

Radiative Transfer Approximations 4 

Asymptotic Theory 4 

Two-stream Approximations 7 

Results 11 

Concluding Remarks 39 

References 40 


v 




COMPARATIVE ACCURACY OF THE ALBEDO, TRANSMISSION AND ABSORPTION FOR 
SELECTED RADIATIVE TRANSFER APPROXIMATIONS 


Michael D. King and Harshvardhan 1 


Laboratory for Atmospheres 
Goddard Space Flight Center, NASA 
Greenbelt, MD 20771 


SUMMARY 

In the present study the plane albedo, total transmission and fractional absorption predicted 
by eight different radiative transfer approximations are compared with doubling computations as a 
function of optical thickness, solar zenith angle and single scattering albedo. The phase function 
used in these computations is intended to be representative of a fair weather cumulus cloud having 
an asymmetry factor g = 0.843. These results show that specific regions can be identified where 
one radiative transfer approximation is more accurate than another. For climate model applica- 
tions in which a single radiative transfer approximation is required to accurately model reflection, 
transmission and absorption at all optical depths, solar zenith angles and single scattering albedos, 
our results show that this requirement is not satisfied by any of the approximate methods. How- 
ever, for remote sensing applications involving flux measurements of either reflected or transmit- 
ted radiation, it is generally possible to use these results as a guide in selecting the most accurate 
approximation to use. 

Of the models we have examined, asymptotic theory (Sobolev, 1975, van de Hulst, 1980), al- 
though seldom used, is the most accurate approximation for optical thicknesses greater than about 
6, where errors less than 5% are obtained for reflection, transmission and absorption for all solar 
zenith angles and all single scattering albedos. Among the two-stream methods, Coakley-Chylek’s 
model 1 is accurate to within 5% for optically thin atmospheres having optical thicknesses less 
than about 0.2 for most values of the solar zenith angle. Though the accuracies of the delta-Ed- 
dington and Meador-Weaver approximations are less easily summarized, it can generally be con- 
cluded that the delta-Eddington approximation is the most accurate two-stream method for con- 
servative scattering in optically thick atmospheres, whereas the Meador-Weaver approximation is 
the most accurate two-stream method for nonconservative scattering (co Q < 0.9). In general, the 
PIFM and delta-discrete ordinates methods are quite similar to, but slightly worse than, the delta- 
Eddington method over the entire range of variables. The other two-stream models presented in 
this report are the Eddington approximation and Coakley-Chylek’s model 2. In both of these 
methods it is always possible to find one of the other six models for which more accurate esti- 
mates of the plane albedo, total transmission or fractional absorption are obtained. This precludes 
the use of either of these models in any radiative transfer problem involving large asymmetry fac- 
tors, such as those examined in this report. 


INTRODUCTION 

In many remote sensing and climate model applications, it is necessary to calculate the plane 
albedo, total transmission and fractional absorption of atmospheric layers as a function of opti- 
cal thickness and solar zenith angle. Among the simplest and most widely used approximations 
to the radiative transfer equation are the two-stream approximations. These approximations 

^ Also affiliated with Department of Meteorology, University of Maryland, College Park, MD 20742. 
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have been discussed and analyzed by Irvine (1968), Kawata and Irvine (1970), Shettle and Wein- 
man (1970), Liou (1973, 1974), Coakley and Chylek (1975), Joseph et al. (1976) and Schaller 
(1979). Recently, Meador and Weaver (1980) and Zdunkowski et al. (1980) have shown that a 
whole class of approximate two-stream solutions can be reduced to a standard form with only a 
few coefficients. These coefficients depend on the cosine of the solar zenith angle (^q), single 
scattering albedo (co Q )and one or more moments of the single scattering phase function, while the 
general equations- for the plane albedo and total transmission depend in addition on the total opti- 
cal thickness of the layer (Tj). In spite of this long history of development, no generally agreed- 
upon variable ranges exist within which one can use a given approximation with assurance that 
accurate reflected, transmitted and absorbed flux densities will be obtained. 

Once a layer is sufficiently thick, an asymptotic regime is established within the layer which 
permits the extension of the plane albedo, total transmission and fractional absorption to com- 
parable values for a semi-infinite layer (van de Hulst, 1968a, 1980; Sobolev, 1975). In the present 
study, this asymptotic method for thick layers is compared in accuracy with selected two-stream 
approximations forO<^ 0 <l and 0.1 <r t <100. Results have been obtained for four values of 
the single scattering albedo (viz., 1.0, 0.99, 0.9 and 0.8) and for a cloud phase function having an 
asymmetry factor g = 0.843. Following the suggestion ofWiscombe and Joseph (1977), who con- 
sideied the accuracy of the Eddington approximation for g^0.5, both absolute and relative errors 
in the plane albedo, total transmission and fractional absorption aie presented in this report. We 
concentrate our comparisons on asymptotic theory for thick layers and Coakley-Chylek’s models 1 
and 2, Meador-Weaver, Eddington, delta-Eddington, PIFM and delta-discrete ordinates approxima- 
tions. 


MULTIPLE SCATTERING COMPUTATIONS 

To provide a baseline for assessing the accuracy of various radiative transfer approximations, numerical compu- 
tations were performed for a model atmosphere composed of cloud particles. Fig. 1 illustrates the phase function 
employed in these calculations, which is based on Mie theory for a wavelength X = 0.754 /am, refractive index 
m = 1.332. and a size distribution of particles of a given radius proportional to r 6 exp(-1.6187r), where r is the 
particle radius in /mi. This distribution of particles is a gamma distribution with an effective radius of 5.56 /am and 
an effective variance of 0.1 1 1. and is considered typical of fair weather cumulus (FWC) clouds (Hansen, 1971), 
This distribution is similar to Deinnendjian’s (1963) cloud C. 1 model, except that the effective radius in Deirmend- 
jian’s model is 6.0 /am. 

In performing our radiative transfer calculations, we have followed the common practice of expressing the prod- 
uct of the single scattering albedo co Q and phase function <I>(cos 0), as a finite expansion of Legendre polynomials of 
the form 


L 

w 0 <P(cosQ) = ^ LO l p^cos 0), (!) 

/= 0 

where 0 is the scattering angle and P } (cos 0) a Legendre polynomial of order /. With this definition, the phase 
function obeys the normalization condition 


1 

'A f <I>(c<as &)d(cos 0) = 1, 

-1 

with the asymmetry factor g related to the Legendre coefficient go. by 


( 2 ) 


/if <l>(cosQ) cosQ d(cosQ) 

-1 • 


8 = /(3 Uq ) • 


(3) 





PHASE FUNCTION 


FWC MODEL 

X = 0.754 /xm 
m = 1 .332 — O.Oi 
a = 0.84276 
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Phase function as a function of scattering angle for a fair weather cumulus (FWC) size distribution given 
by n(r) exp{— 1.6187r), where X = 0.754 pm and m = 1.332 - O.Oi. 




The asymmetry factor for the FWC cloud model is g= 0.843. 

Multiple scattering calculations were performed for the azimuth-independent term of the reflection and trans- 
mission functions using the doubling method described by Hansen and Travis (1974), together with the invariant 
imbedding initialization described by King (1983). In terms of these functions, the azimuthally-averaged reflected 
1 (0, -/r) and transmitted 7 u (r t , y) intensities from a horizontally homogeneous atmosphere illuminated from above 
by a parallel beam of radiation of incident flux density F Q may be expressed as 


/°( 0 ,-/ 0 - (p 0 FJtt)r°(t v ii , p fl ). 

( 4 ) 

7 Vp F) = (P 0 Fq/tt) r°(r t 

( 5 ) 


where 7, is the total optical thickness of the atmosphere, the cosine of the solar zenith angle and /a the cosine of 
the zenith angle with respect to the outward normal (0 <ai, F 0 < 1). 

In terms of the azimuth-independent reflection function R°(r t ;ju, n ) and transmission function T°( r. ;/i,p n ), 
the plane albedo [Kz" t , /i Q )] , total transmission [r(7 { . /!„)] and fractional absorption [a(r t . p )] of the layer arc 
given by 10 

1 

r( T v V 0 ) = 2 J R°(T t \n,n 0 )ndn, ( 6 ) 

0 

1 

= 2 j7’°(T t .7i,p 0 )^r7/i + exp(-7 t //a 0 ), (7) 

0 


a ( 7 t ’^o) “ 1 - r (T t ,V 0 ) -t(T v n Q ) . (8) 

Due to the use of a highly anisotropic phase function and the requirement that accurate computations of the plane 
albedo, total transmission and fractional absorption be obtained, we have subdivided the u angular interval using 
a Gaussian quadrature of order 80 (King, 1983). 

Since the major purpose of this study is to examine the accuracy of various radiative transfer approximations 
over a wide range oi optical depths, solar zenith angles and single scattering albedos, we have ignored the effects of 
surface reflection. As a result, the reflection and transmission functions appearing in (4) - (7) apply to those of an 
isolated cloud layer only. 


RADIATIVE TRANSFER APPROXIMATIONS 


Asymptotic Theory 

When the optical thickness is sufficiently large, an asymptotic regime is established within the layer such that 
the reflection and transmission functions can be expressed in terms of functions applicable to a semi-infinite layer 
(van de Hulst. 1968a, 1980). From these expressions, coupled with the definitions of plane albedo and total trans- 
mission given previously, it can be shown that the asymptotic theory approximations for the plane albedo [r(r ,ju )] , 
total transmission [r(Tj. Pq)] and fractional absorption [ 0 ( 7 ^. Mq)] of the layer are given by 


?( 7 t . v 0 ) = rjfi 0 ) - 


m n l 


l-/2 e -^t 


K(n Q )<r 2kT t 


( 9 ) 
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m n 


( 10 ) 


?(T t ,H 0 ) = 


l-/2e-2*t 


A'(m 0 ) 


„-kr. 


&(r t ,P 0 ) = 1 -r(r t ,/u 0 )-?(T t ,n 0 ) 


(ID 


In these expressions r^f j q ) is the plane albedo of a semi-infinite atmosphere, A(y Q ) the escape function, k the 
diffusion exponent describing the attenuation of radiation in the diffusion domain, and m, n and / constants which 
depend primarily on the single scattering albedo and asymmetry factor (King, 1981). 

The escape function and diffusion exponent, as well as other asymptotic functions and constants appearing in 
(9) and (10), can be obtained by applying the asymptotic fitting method of van de Hulst (1968b). In this method, 
computational results from the doubling method are fit to asymptotic expressions for the plane albedo, diffuse 
transmission and internal intensity field as a function of the optical depth. The functions K(P 0 ) and thus 

obtained are illustrated in King and Harshvardhan (1986) as a function of for the FWC phase function and for 
four values of co 0 . In addition to these functions, the plane albedo and total transmission of thick layers depend on 
the constants m, n, l and k [cf. Eqs. (9) and (10)] . Each of these constants is strongly dependent with a some- 
what weaker dependence on g. King (1981) has shown that m, n and / can be well described by a function of a sim- 
ilarity parameters, defined by 



( 12 ) 


where s reduces to (1 - co 0 )^ for isotropic scattering and spans the range 0 (cj q - 1) to l(w 0 - 0). 

Although the diffusion exponent k does not obey such a similarity relationship, the function fc/( 1 - w 0 g) does. 
Fig. 2 illustrates kj( 1 - co Q g) as a function of s for both the FWC (co 0 - 0.99, 0.9 and 0.8) and Henyey-Greenstein 
phase functions for varying values of (0.9999, 0.999, 0.996, 0.99, 0.96, 0.9, 0.8 and 0.6) and g (0.8, 0.85 and 
0.9). The computational results presented in Fig. 2 were fit to k/ (1 — co Q g) as a function of s. The formulas for 
m, n, / and kj{ 1 — w^g) are summarized in Table 1. The formula for m is identical to that obtained by King (1981), 
whereas the coefficients in the formulas for n and l differ slightly in order to give a better fit for small values of s. 

For conservative scattering, when r^/^) = n = / = 1 and m = k = 0, the asymptotic expressions for the plane 
albedo and total transmission given in (9) and (10) are indeterminate. Expanding n, l, m and k to first order in s. it 
can be shown that (9) and (10) can be rewritten as (King, 1981) 


r{r v H Q ) 


f ( T f 


4A(p 0 ) 

■ - - > 

3(1 -g)(r t + 2q 0 ) 

4 K(p 0 ) 

3(1 -g)(T t + 2 q 0 ) 


(13) 

(14) 


where is the extrapolation length. The reduced extrapolation length q = (1 - g)q Q is known to range between 
0.709 and 0.715 for all possible phase functions (van de Hulst , 1980), and can be well approximated by 0.714 for 
anisotropic cloud phase functions (King, 1981). Thus it is seen that the plane albedo and total transmission in opti- 
cally thick, conservatively scattering layers are a function of ( 1 -g)r t , with all of the solar zenith angle dependence 
contained in K(fJ Q ). 
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HENYEY-GREENSTEIN 
o g = 0.80 
• g = 0.85 
□ g = 0.90 

FWC 

■ g = 0.84276 


Figure 2. kl( 1 as a function of similarity parameter s= [(1 - cj q )/(1 -co Q ,g-] ” , where k is the diffusion 

exponent. The symbols represent values obtained by numerical computation for FWC and Henyey- 
Greenstein phase functions, and the curve the result of a least-squares fit to an analytic equation. 






TABLE 1 


Similarity relations satisfied by the constants which arise in asymptotic expressions for the plane albedo, 
total transmission and fractional absorption of thick layers. 


(1 —0.681 s)(l -s) 
(1 + 0.792s) 


n 


(1 + 0.414s) (1 -s) 
(1 + 1.888s) 


a 


m = (1 + 1.537s) In 


1 + 1.800s- 7.087s 2 + 4.740s 3 
(1 — 0.819s)(l -s) 2 


*/(l = s - 


(0.985 - 0.253s)s 2 
(6.464 - 5.464s) 


Two-Stream Approximations 

The two-stream approximations in radiative transfer are based on assuming various analytic forms for the up- 
ward and downward intensity fields within and at the boundaries of a plane-parallel medium. Substituting the as- 
sumed angular distribution into the integro-differential form of the equation of transfer results in a set of differen- 
tial equations for the upward [F~(t. f 0 ) 1 and downward [F + (t, /i Q )] flux densities (Meador and Weaver, 1980; 
Zdunkowski et al, 1980); 

dF (T :^ = 

(It 

= 7 2 F (r, M 0 ) -7|^ + (FFo) +F 0 cj 0 t 4 e , 

dr 


(15) 

( 16 ) 


where 


1 

2n f /°(r,±M)M^. 
0 


(H) 


In order to obtain the forms given in (15) and (16), it is often necessary to approximate the scattering phase 
function in order to integrate the source function analytically. Expressions for the plane albedo, total transmission 
and fractional absorption are obtained by solving (15) and (16). subject to the boundary conditionsF (r t . y 0 ) 
F + (0, ) = 0. The results may be obtained in the form (Meador and Weaver. 1980) 
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r(r v ^ Q ) 


coo 

(1 -A^jHOc+T,)^ +(k- 7,)e-' Cr t] 


x [(1 -fcM 0 )(o 2 +/C73)e kT t 


-(1 +fcM 0 )(a 2 “fc7 3 )e kT t-2k(y 3 -a 2 n Q )e~ T t/^O], 


( 18 ) 


?(r t ,/i 0 ) = e r t/^ 0 x | 1 — 

( (1 [(/c + 7 j)e^ T t + (fc — 7])e~^ Tt ] 


x [(1 +fcp 0 )(a 1 +/c7 4 )e fcTt -(l - kn Q )(a ( - ky^)e~ kT i 

-2A:(7 4 +«! M 0 )e r t^o] j , (19) 

where 

a l =7,7 4 +7 2 7 3 , (20) 

a 2 = 7iT 3 +7 2 7 4 , (21) 

K wi 7 2 ) , (22) 

y 4 =1 ~ y 3' (23) 


and the fractional absorption a(r t ,/a 0 ) is given by (1 1). 

The 7 1 ■ 7 2 and 7 3 coefficients in (1 5) and (16) for various two-stream approximations, along with references to 
their original description in the literature, are given in Table 2. Several of these methods employ delta scaling (Jo- 
seph et al.. 1976) in which a fraction / ot the scattered energy is considered to be in the forward peak, approximated 
as a Dirac delta function. For each of these methods, which include the delta-Eddington, Practical Improved Flux 
Method (PIFM) and delta-discrete ordinates methods, (18) and (19) can still be used as long as the following trans- 
lormations are made in the coefficients and solutions: 

= (1 -w 0 i)r t , (24) 


w O _>w o' = ( 1 ~/) w o/0 (25) 

g~*g’ = (g- Ml ~f). (26) 


The primed quantities in (24) - (26), when substituted into the expressions for 7 j , y 2 and 7 3 (cf. Table 2) as well as 
into (18) and (19), yield the relevant expressions for r(T { , n Q ), ?(r,, M 0 ) and a(r v p Q )for the delta-scaled approxima- 
tions. Though various choices of / are possible, the most frequently used choice, and the one used in all computa- 
tional results to be presented below, is f=g 2 . 
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TABLE 2. Summary of 7 . coefficients in selected two-stream approximations. 


Method 

Reference 

*1 

*2 


Eddington 

Kawata and Irvine (1970) 

1/4(7 -co () (4 + 3s)] 

-1/4(1 -co 0 (4 -3s)] 

1/4(2 -3gu Q ) 

delta-Eddington 

Joseph et al. (1976) 

1/4(7- gj (1 '(4 + 3s')] 

-1/4(1 -w 0 '(4 - 3s0] 

1/4(2 - 3 g’n Q ) 

PIFM 

Zdunkowski et al. (1980) 

1/4(8 -oj 0 '(5 +3s')l 

3/4(co 0 '(l -/)! 

1 / 4(2 - 3s'/i 0 ) 

discrete ordinates 

Liou ( 1973, 1974) 

■\JT / 2 [2 - or 0 (l + s)l 

>/3/2[w 0 (l -S)] 

1/2(1-^3sm 0 ) 

delta-discrete ordinates 

Schaller (1979) 

^/2[2 - w 0 '(l +^ 0 ] 

yjJ / 2 [u q '(\ -s')] 

1 / 2(1 - g'ny) 

Coakley-Chylek (1) 

Coakley and Chylek (1975) 

{1 - w () [l 

u oP(» 0 )l» 0 

P<M 0 ) 

Coakley-Chylek (11) 

Coakley and Chylek (1975) 

2(1 - w Q (l - 0 )] 

2 utf 

P ( V 0 ) 

Meador-Weaver 

Meador and Weaver ( 1 980) 

7 -3g 2 - w Q (4+3g)+ w 0 g 2 [4/3(p 0 )+3 g] 

- 1 3 g 2 fw 0 (4 - 3s) + ^()g 2 [40(a< o )+3s-4] 


4(1 -S 2 (l — M 0 )l 

4(1 -s 2 (l-P 0 )] 



fined a S the Mead ° r ' WeaVer and Coakl ey-Chylek approximations, 7 3 is set equal to the backscatter fraction 0^). de- 


/5(h 0 ) - f h°(-n,u Q )dn, 

2co J 


(27) 


’o o 


where h°( M v Q ) is the azimuthal average of co 0 <D( cox 0) for incident solar radiation in the direction ju n and re- 
fiected m the direction - u. The backscatter fraction 0 (m o ) is illustrated in King and Harshvardhan (1986) for the 
rWC phase function used in the present investigation. 

The Coakley and Chylek (1975) model 2, which we will hereafter refer to as Coakley-Chylek (II) uses the 
average backscatter fraction 0 in the expression for and 7 2 . This constant is defined as 


P= /^ 0 ) rf % 

0 

For the phase function used in the present investigation,^ = 0.1772. 

For conservative scattering, for which y x = J 2 and k = 0. Eqs. (18) and (19) reduce to 

>V t ,M 0 )= 1 -?(T t ,p 0 ) = — f7i r t + (7 3 -7 1 // n )(I -e~ T t l^o)]. 

7 V t 


(28) 


(29) 


For Coakley and Cheek's (1975) model 1 . hereafter referred to as Coakley-Chylek (I), 7 3 =7, p n and thus P9) re- 
duces to an especially simple form (see Table 2). ' M 0 Jna tnus re 

Since k exceeds unity for strongly absorbing atmospheres in all two-stream approximations (cf. Table ■>) con- 
itions can easily exist tor which u 0 = k ~\ especially in the water vapor bands. Though this condition can lead 

a numencal singularity m (18) and (19) the singularity is removable, and when » = *“ 1 it is rather straight- 
forward to show that (18) and (19) reduce to 0 straignt 

co 0 


r(T V^ 


- X j( a 2 M 0 +7 3 )e Tt/ ^0 - [(a 2 /J 0 + 7 3 ) 

[(1 + 7 , M 0 ) c7t ^0 +(1 - 7 ] M 0 ) e ' Tt ^ 0l 

+ 2<a 2 „ 0 - T3 )r t /„ 0 l e-n /« 0 | 


(30) 


r(r t ,p () ) = e T t/^o x 


1 - 




2 [( 1 +7j V 0 )e T t/t A c , + 0 -7,,h 0 )e f t/ p 0] 
x |K a i M 0 -7 4 )-2(a, Mo + T 4 )7 t /p 0 ]e 7 t/M ( , -(a. M fl ~ 7 4 ) e~ r t/M 0 ] j 


(31) 


This 


use may be avoided by either applying these formulae when n = k~\ or by displacing p n by a very small 
increment and applying (18) and (19). as suggested by Zdunkowskief al. (1980). 

In computing fluxes for multi-layer systems overlying a reflecting surface, it is also necessary to compute the 
albedo and transmission of layers for diffuse radiation. For parameterization purposes, it is usual to compute these 
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quantities for an isotropic incident source. Under this situation, the global (spherical) albedo and global trans- 
mission can be obtained by integrating the corresponding plane albedo and total transmission solutions as a function 
of n Q : 1 ^ 

7(r t ) = 2 / r(T t ,H 0 )n 0 dn 0 , 

0 

1 

T (r t ) = 2 J?(T t ,*i 0 ) - 

0 

Inspection of (18) and (19), or even the simpler (29) for conservative scattering, shows that a general closed form 
solution for (32) and (33) does not exist. Although these integrations can be carried out numerically for each 
specific two-stream model, this is not practical for most modeling applications. For those two-stream models for 
which 0(n ) does not explicitly appear in any of the coefficients, it is in principle possible to analytically integrate 
(18) and (19) to obtain expressions for r(Tf) andT^lfor specific models. In asymptotic theory it is possible to ob- 
tain simple analytic formulae for F(r t )and f(T t ). These expressions can be found in King (1981). 

Coakley and Chylek (1975) suggest that (32) and (33) may be avoided by using a second set of globally-aver- 
aged two-stream equations in which the incident isotropic radiation is treated as an upper boundary condition. In 
general, however, the results obtained by this approach yield different results from those obtained using (32) and 
(33) for the same model. In the results to be presented below, we have restricted our comparisons of the radiative 
properties of various radiative transfer models to those for an isolated layer composed of cloud particles. The com- 
plicating effects of ground albedo have therefore not been considered, but a comparable intercomparison of results 
for r(T t ) and 7(r t )as a function of r t and co 0 remains an important study to be pursued in the future. 


(32) 

(33) 


RESULTS 

We have examined both the absolute and relative accuracies of the plane albedo, total transmission and frac- 
tional absorption (where applicable) as a function of Tj and Hq for four values ofthe single scattering albedo (cOq - 1 .0, 
0.99. 0.9 and 0.8) and for all radiative transfer approximations discussed in the previous section. Fig. 3 illustrates 
r(r . H 0 ) = 1 - t{r v H 0 ) as a function of r t (0.1 <r t <100) and H 0 (0 <u Q < 1) for conservative scattering (co 0 = 1 ), 
and Fig. 4 illustrates corresponding results for r(r t , H 0 ) for each of eight different radiative transfer approximations. 
The shaded regions in the asymptotic theory and Eddington approximation results of Fig. 4 delineate regions for 
which nonphysical (negative) plane albedos are obtained in the approximations. In both of these models, nonphysi- 
cal regions occur for optically thin atmospheres and for small solar zenith angles (large values of H 0 )- 

Fig. 5 illustrates a 4 x 3 plot composite of errors for conservative scattering and for four radiative transfer 
models, where the first row applies to asymptotic theory and succeeding rows to thedelta-Eddington, Meador-Weaver 
and Coakley-Chylek (1) approximations. Individual plots in the first column of Fig. 5 represent absolute errors in 
the plane albedo, defined as 


Az(t , ,M 0 ) = r(T v n Q )-r(T v fi Q ). 


(34) 


With this definition, positive (negative) errors indicate that the radiative transfer approximation overestimates 
(underestimates) the exact albedo, taken as the computational results presented in Fig. 3. Similar definitions apply 
to errors in the total transmission [Ar(r t , H 0 )] and fractional absorption [Aa(T t , h 0 )]- The relative errors in the 
plane albedo [Ar(r r H 0 )Mt v h 0 )] and total transmission [Ar(r t , H 0 )/f(r t ,H 0 )] are presented in succeeding col- 
umns of Fig. 5, and are given in per cent. Relative errors with magnitudes greater than 20% and absolute errors 
with magnitudes greater than 0.2 are not plotted in Fig. 5 or in subsequent figures. 
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ASYMPTOTIC THEORY EDDINGTON 

T(t,)1 0 I T~It.Ho) OJQ = 1.0 




Ho Mo 


Figure 4. Computations of the plane albedo as a function of optical thickness and cosine of the solar zenith 
angle for eight different radiative transfer approximations identified in the figure. 
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Figure 5. Absolute and relative accuracy of asymptotic theory, delta-Eddington, Meador-Weaver and Coakley- 
Chylek (I) approxmations to the plane albedo and total transmission as a function of optical thickness 
and cosine of the solar zenith angle for conservative scattering (w Q = 1 . 0 ). All relative accuracy values 
are in per cent. The FWC phase function is assumed throughout. 
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Fig. 5 shows that asymptotic theory is accurate to within 5% in both reflection and transmission for ~>3 

when <0.9. Furthermore, asymptotic theory is accurate to within 1% for all solar zenith angles when >»8. 

The Coakley-Chylek (I) approximation, on the other hand, is accurate to within 5% for all < 0.2 when n 0 >0.\. 
Since r(r t , n 0 ) is small when r t /M 0 < 0.5 (cf. Fig. 3), a relative error of 5% is too stringent a criterion to use for 
accepting a model in this range. Adopting instead the absolute error criterion I Ar(r t , M 0 ) I < 0.005, we see 
that the range of validity of the Coakley-Chylek (I) approximation can be extended to r t <0.5 for a wide range 
of solar zenith angles. An advantage of both of these models is that their absolute and relative accuracies show 
little sensitivity to /i Q in their respective ranges of validity. 

In contrast to these models, the delta-Eddington approximation tends to have its greatest accuracy when 
H Q > 0.5, regardless of optical thickness. The large relative albedo errors which occur when r t //i Q is small are not 
critical, since the absolute errors are small in this range. Similarly, when t ^ 1> 10 and <0.5, the small values of 
f(r t , ju Q ) allow one to extend the range of validity of the delta-Eddington model to values of n 0 somewhat lower 
than 0.5. 

The Meador-Weaver approximation, which was developed as a composite of the Eddington and Coakley-Chylek 
(I) methods, has most of the characteristics of the latter for conservative scattering, especially for thick atmo- 
spheres. In optically thin atmospheres, on the other hand, it has a much greater sensitivity than the Coakley- 
Chylek (I) method. This makes the Meador-Weaver approximation less suitable than alternative methods for con- 
servative scattering over the entire range of variables, at least for the high values of asymmetry factor considered 
in the present investigation. 

Detailed results for conservative scattering analogous to Fig. 5 are presented in Fig. 6 for four other radiative 
transfer approximations, viz., the Eddington, Coakley-Chylek (II), PIFM and delta-discrete ordinates methods. On 
examination of Table 2, one can readily show that the PIFM method of Zdunkowski ct al. (1980) is identical to 
the delta-Eddington method for conservative scattering, which accounts for the identical appearance of the corre- 
sponding panels in Figs. 5 and 6 for these two models. The difference Zdunkowski et al. report between the delta- 
Eddington and PIFM methods when = 1 is likely a result of their using different values of / in the scaling for- 
mulae (24) - (26) for each method. 

A careful examination of the results in Fig. 6 shows that the Eddington, PIFM and delta-discrete ordinates 
methods are quite accurate for optically thick, conservative atmospheres, although they differ substantially in 
their accuracy at small and intermediate optical depths. From a comparison of the results in Figs. 5 and 6, it fol- 
lows that the Coakley-Chylek (II) approximation is everywhere less accurate than the Coakley-Chylek (I) approx- 
imation. In general, it can be concluded that when co Q = 1 the majority of models considered in this report are in- 
accurate in their approximation of the plane albedo when jUq < 1 and 1 <r^. <5, with asymptotic theory being 
the best suited in this difficult range of variables. 

Fig. 7 illustrates doubling computations of the plane albedo, total transmission and fractional absorption as a 
function of t. and / 1 q for nonconservative scattering (Wq = 0.99), where we have used the same phase function il- 
lustrated in Fig. 1, but simply scaled the Legendre coefficients by Wq. Figs. 8 and 9 show corresponding figures for 
each of the eight radiative transfer approximations considered in this report. Again we have used shading to iden- 
tify all nonphysical regions where cither the plane albedo or fractional absorption is negative or the total transmis- 
sion exceeds unity. 

Fig. 10 illustrates a 4 x 3 plot composite showing absolute errors in the plane albedo, total transmission and 
fractional absorption as a function of and / 1 q for Uq = 0.99 and for each of the four models presented in Fig. 8. 
A similar plot composite for the remaining four models is presented in Fig. 11. Corresponding results for relative 
errors for all eight models are presented in Figs. 12 and 13. In all of these figures, individual plots in the first col- 
umn represent errors in the plane albedo, while plots in succeeding columns represent errors in the total transmis- 
sion and fractional absorption, respectively. Due to the small amount of absorption when Wq - 0.99 and < 5, 
the relative and absolute errors for all eight models are nearly identical to corresponding errors in the conserva- 
tive case. Only in the case of high sun and large optical depths do appreciable changes occur in the accuracy of dif- 
ferent models. Further discussion of these differences is reserved for the following treatment of the cOq - 0.9 case. 

Fig. 14 illustrates doubling computations of the plane albedo, total transmission and fractional absorption for 
cj q = 0.9, while Figs. 15 and 16 show corresponding figures for all eight radiative transfer approximations consid- 
ered in this report. Figs. 17-20 illustrate absolute and relative errors obtained for co Q = 0.9 in the same format 
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as used previously in Figs. 10 - 13. Figs. 17 and 19 show that asymptotic theory is equally as valid an approxi- 
mation for optically thick, nonconservative atmospheres as it is for optically thick, conservative atmospheres (cf. 
Fig. 5). Relative errors of 5% or less are achieved in asymptotic theory for reflection, transmission and absorp- 
tion when x t >6, regardless of solar zenith angle. For cases in which reflection is the most important, the results 
presented in Fig. 19 show that asymptotic theory can be applied to optical depths as low as 4 with an accuracy 
of better than 5%. 

As in the case of conservative scattering, the Coakley-Chylek (I) method is the most accurate approxima- 
tion for optically thin atmospheres, with a tendency to be somewhat more accurate for small solar zenith angles 
(large values of ju Q ). In order to have an accuracy of better than 5% in reflection, transmission and absorption, 
Fig. 19 suggests that it is necessary for r t <0.1 and n Q >0.1. However, since a( 7 t ,p 0 )< 1 when x t //ig <0.5 (cf. 
Fig. 14c), it is more meaningful to use the absolute error criterion |Aa(r t , /ig)| < 0.005. Thus we conclude that 
the range of validity of the Coakley-Chylek (I) approximation can generally be extended to include all optical depths 
less than some maximum in the range 0.2 <x t <0.7, depending on solar zenith angle. 

For the delta- Eddington approximation, comparison of Figs. 5 and 19 show that relative errors in the plane al- 
bedo and total transmission degrade somewhat in accuracy as absorption increases, especially for optically thick atmo- 
spheres. This is a consequence of the fact that the ^(Mq) computed in the delta-Eddington approximation is nearly 
linear in for all single scattering albedos, whereas the true /^(jUg) has increasing curvature as absorption increases 
(cf. King and Harshvardhan, 1986). Both the absolute and relative errors in the other delta-scaled approximations 
(viz., the PIFM and delta-discrete ordinates methods) are very similar to, but slightly worse than, the delta-Eddington 
method. A feature of all of these methods is the isolated region at intermediate values of x t and p 0 where absorp- 
tion errors in excess of 10% occur (cf. Figs. 19 and 20). 

Although the Meador-Weaver approximation was previously shown to be an inferior model for conservative scat- 
tering, it is clear from Figs. 17 and 19 that its accuracy improves dramatically as absorption increases, especially for 
reflection. This is true for both optically thin and thick atmospheres. Moreover, it is the only two-stream model 
which has an albedo accuracy of better than 5% over a wide range of solar zenith angles when r t >2, although the 
absorption error is sometimes as large as 10% in this range of variables. The explanation for the exceptional accur- 
acy of the Meador-Weaver approximation in optically thick, nonconservative atmospheres is that exhibits 

significant curvature in for all single scattering albedos, as does the true (ju Q ) (cf. King and Harshvardhan, 
1986). This feature is unique to the Meador-Weaver method among two-stream approximations. 

Although the Eddington approximation was previously shown to be quite accurate for optically thick, nearly 
conservative atmospheres, the results presented in Figs. 18 and 20 show that the Eddington approximation has no 
useful regime where albedo errors less than 5% occur when'. co 0 = 0.9. From a similar examination of results for the 
Coakley-Chylek (II) method, we conclude that this model is also a poor approximation for nonconservative atmo- 
spheres involving collimated radiation. 

Detailed results analogous to Figs. 14 - 20 are presented in Figs. 21-27 for to 0 = 0.8. From a comparison of 
the relative errors in Figs. 5, 12, 19 and 26, one can see that there is a clear tendency for the Meador-Weaver approx- 
imation to improve in accuracy as cjg decreases, at least to = 0.8, at which point albedo errors of less than 7.5% 
occur for all optical depths when > 0.2. In addition, these figures show that asymptotic theory continues to be 
a valid approximation for optically thick, nonconservative atmospheres. Although the Coakley-Chylek (I) model can 
still be used for optically thin atmospheres, the other two-stream methods generally degrade in accuracy as absorption 
increases. Again it is clear from Figs. 26 and 27 that all of the delta-scaled two-stream methods have an isolated 
region at intermediate values of x t and ju Q where absorption errors in excess of 10% occur. In general, it can be con- 
cluded that all models with the exception of asymptotic theory and the Meador-Weaver approximation have diffi- 
culty in optically thick, nonconservative atmospheres. 
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Figure 21. Doubling computations of the (a) plane albedo, (b) total transmission and (c) fractional absorption as a 
function of optical thickness and cosine of the solar zenith angle for a FWC phase function with noncon 
servative scattering (co Q = 0.8). 
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CONCLUDING REMARKS 


Most previous intercomparisons of radiative transfer approximations have concentrated on presenting results 
for the plane albedo as a function of cosine of the solar zenith angle for selected values of the optical depth. On 
some occasions, intercomparisons have been further restricted to selected values of Although plane albedo er- 
rors in the delta-Eddington approximation are less than 5% for all optical depths when Wq = 0.9 and = 0.4 
(cf. Fig. 19), a generalized conclusion on its overall accuracy based on this restricted intercomparison would be 
highly misleading. We therefore feel it is important to examine the accuracy of the plane albedo, total transmis- 
sion and fractional absorption as a function of optical depth and solar zenith angle before drawing conclusions 
about the overall accuracy of a given approximation. 

After examining a wide variety of radiative transfer approximations over a large range of optical depths, solar 
zenith angles and single scattering albedos, it has become evident why some approximations succeed while others 
fail in specific regimes. For example, a straight-forward comparison of (13) and (29) shows that when co Q = 1 and 
r t /n 0 ^ 1, the plane albedo obtained from asymptotic theory and two-stream approximations are equivalent, pro- 
vided the two-stream coefficients y^ and ? 3 satisfy the following criteria: 

*1 = “(1 -*). (35) 

2 K(/i n ) 

*3 " 1+ ^o ~ • (36) 

3 q' 

For the delta-scaled approximations, the only difference in these criteria is the substitution g -*■ g' in (35). Since 
q ~ 0.714 for all possible phase functions, (35) implies that ~ 0.7 (1 — g ) for unsealed approximations and 
0.7 (1 — g*) for scaled approximations. Table 2 shows that the Eddington and delta-Eddington methods satisfy these 
requirements the most closely. Furthermore, in order for a two-stream method to perform well for optically thick, 
conservative atmospheres, it is necessary for y^ to be a linear function of over most of the range of solar zenith 
angles. This readily follows from a comparison of (36) with the escape function K(p 0 )presented in King and Harsh- 
vardhan (1986, Fig. 4). The poor performance of the Meador-Weaver and Coakley-Chylek (I) methods under these 
conditions is at least in part a result of their choice of y 3 = /3 (ju q ), a function which is highly nonlinear in /r Q . The 
Eddington and delta-Eddington methods, on the other hand, very nearly satisfy (36). Note that although 7 3 can 
be negative for high sun in the Eddington method, this is in accord with the requirement given by (36). 

Due to the importance of developing a radiative transfer approximation which is accurate for all solar zenith 
angles and over a wide range of optical depths, our results suggest that a hybrid two-stream model that reduces to 
asymptotic theory for thick atmospheres but extends the range of validity of asymptotic theory to thinner atmo- 
spheres would be extremely valuable. The development of such a model remains a challenge for further study. 
Finally, we would like to note that none of the conclusions drawn in the present investigation are affected by our 
choice of a Mie theory phase function. Limited intercomparisons with doubling computations using the Henyey- 
Greenstein phase function with the same asymmetry factor as in the fair weather cumulus model ( g = 0.843) yield 
error plots with virtually the same appearance as those of the FWC phase function. 


Goddard Space Flight Center 
National Aeronautics and Space Administration 
Greenbelt, MD 20771 
September 6, 1985 
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